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ABSTRACT 

We discuss the origin of the Tully-Fisher (TF) relation using the 7V-body/SPH method, which includes 
cooling, star formation and stellar feedback of energy, mass and metals. We consider initially rotating 
overdense spheres, and trace formation processes of disk galaxies from z — 25 to z = in the Cold 
Dark Matter (CDM) cosmology. To clarify the origin of the TF relation, we simulate formation of 14 
galaxies with different masses and spin parameters, and compute observable values, such as the total 
magnitude and the line-width. We find that the simulated galaxies reproduce the slope and scatter of 
the TF relation: the slope is originated in the difference of total galactic masses, and the scatter is 
produced by the difference of initial spin parameters. As well as the TF relation, observed features of 
spiral galaxies, such as the exponential light-profile and the flat rotation curve, are reproduced in our 
simulations, which were assumed a priori in past semi-analytical approaches. 

Subject headings: galaxies:formation — galaxies:evolution — galaxies: kinematics and dynamics — 
galaxies:statistical 



1. INTRODUCTION 

The Tully-Fisher relation (hereafter TF: Tully & Fisher 
1977) is one of the most basic relations in spiral galax- 
ies, and would be a clue to understand the origin of disk 
galaxies. The TF relation has been roughly considered 
as a product of the virial theorem and a nearly constant 
mass-to-light ratio since its discovery. Recently, the origin 
of the TF relation is discussed by, e.g., Silk (1997) and 
Mo, Mao & White (1998). 

In their semi-analytical approach, Mo, Mao & White 
(1998) succeeded in reproducing the TF relation, assum- 
ing a constant mass-to-light ratio and empirical profiles 
of disks and haloes. Heavens & Jimenez (1999) took the 
same approach but including an empirical star formation 
model, and reproduced the TF relation in four passbands 
simultaneously. In these semi- analytical approaches, how- 
ever, observed features of galaxies, such as the exponential 
profile and flat rotation curve, were not constructed as the 
results of simulations, but assumed a priori. Steinmetz 
& Navarro (1999) performed direct simulations of galaxy 
formation within the cosmological context, and also suc- 
ceeded in getting the slope and scatter of the TF relation. 
They treated the volume much larger than the scale of 
galaxies, and considered environmental effects (e.g., tidal 
field and infall/outflow of mass). However, due to the com- 
plicated behaviors of these effects, it still remains unknown 
what physics produces the TF relation. 

Based on these arguments, we simulate formation and 
evolution of galaxies with the Y-body/SPH method in- 
cluding cooling, star formation and their feedback of en- 
ergy, mass and metals to the ISM. Using the similar 
method, Katz (1992) and Steinmetz & Muller (1994, 1995) 
succeeded in constructing internal structures as observed 
in spiral galaxies, e.g., the exponential density profile, flat 
rotation curve, and distributions of stellar age and metal- 
licity in the bulge, disk and halo. In contrast to their sim- 
ulations of single galaxy formation, we consider formation 



of many galaxies with different masses and spin parame- 
ters. In order to investigate what physics produces the TF 
relation, we concentrate on formation of isolated galaxies 
with different initial conditions. We note that as well as 
the TF relation, observed features of spiral galaxies, such 
as the exponential light-profiles and flat rotation curve, 
are reproduced in our simulations, which were assumed a 
priori in past semi-analytical approaches. 

2. METHODS 
2.1. Numerical Methods 

We use a GRAPE-SPH code, a hybrid scheme of the 
smoothed particle hydrodynamics (SPH) and the TV-body 
integration hardware GRAPE-3 (Sugimoto et al. 1990), 
the code which was firstly discussed by Steinmetz (1996) 
and developed further by us (Koda, Wada & Sofue 1999). 
The code is implemented to treat three-component sys- 
tems: the gas, stars and dark mater, and calculate gravita- 
tional and hydrodynamical forces, which are symmetrized 
to satisfy the Newton's third law. We consider the ra- 
diative and inverse Compton cooling of the gas which is 
optically thin and in collisional ionization equilibrium of 
H and He (X H = 0.76, X He = 0.24; see Katz et al. 1996). 
We also take phenomenological models of star formation 
and their feedback to mimic real galaxy formation. 

The star formation model adopted is basically the same 
as those discussed in Katz (1992), Navarro & White (1993) 
and Steinmetz & Muller (1994, 1995). Stars are formed in 
regions which are locally contracting and Jeans-unstable 
at a rate given by p\ = c*p gas / max(r dyn , r cool ). Here, 
P*i Pgns, Tdyn and t coo i are the densities of stars and the 
gas, local dynamical and cooling timescales, respectively. 
In most cases, Td yn is larger than r coo i, and the star for- 
mation timescale is typically ~ 20rd yn under the adopted 
parameter c* = 0.05. When one third of the mass in a gas 
particle is transformed into stars, we create a new collision- 
less star particle, which inherits the position and velocity 
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of its parent gas particle. 

We take a simple model of feedback. Massive stars with 
> &Mq are assumed to release energy, mass and metals 
into the surrounding gas at a constant rate through stellar 
wind (SW) and type II supernovae (SNII). All of the mas- 
sive stars then become white dwarfs, and 15% of them re- 
sults in type la supernovae (SNIa, Tsujimoto et al. 1995). 
In Table ||, we list total released energy, mass and metal 
per star, and the periods of feedback (see Nomoto et al. 
1997a & 1997b; Yoshii et al. 1996). We can estimate the 
number of stars with > 8M Q , in a star particle, with the 
initial mass function (IMF) of Salpeter (1955). Released 
energy is provided into the surrounding gas as thermal en- 
ergy. Since stars are formed and release energy at dense 
gas regions, where feedback energy is soon radiated away, 
feedback does not so much affect the simulations in our 
model. 

2.2. Initial Conditions 

We simulate formation and evolution of galaxies from 
z = 25 to in the CDM cosmology (fi = I, h = 0.5). 
Our initial conditions are similar to those of Katz (1992) 
and Steinmetz & Miiller (1994, 1995), but we consider 
14 isolated spheres, on which small scale CDM fluctua- 
tions are superimposed with the Zel'dovich approximation 
(Zel'dovich 1970). The densities of the spheres are en- 
hanced above the background field by 5p/p — 0.25. We 
normalize the CDM spectrum so that the rms fluctuation 
in a sphere of radius 8ft. -1 Mpc becomes equal to erg = 0.63 
at z = 0. The spheres are rigidly rotating, and following 
the reduced Hubble expansion (see Steinmetz & Miiller 
1995). Two free parameters, mass and spin parameter, 
are listed in Table @. 

The gas and dark matter are represented by the same 
number of particles, and their mass ratio is set to 1/9. 
The initial temperature of the gas is set to 70 K of the 
cosmic microwave background at z = 25. The mass of 
a gas particle varies between 2.4 x 10 6 and 1.9 x 10 7 M Q 
according to the system mass considered (Table |). The 
mass of a dark matter particle varies between 2.1 x 10 7 
and 1.7 x 10 8 M©. Low resolution may cause an artifi- 
cial heating due to two-body relaxation between the SPH 
and dark matter particles, however this range of particle 
mass is small enough to exclude the artificial heating effect 
(Steinmetz & White 1997). The gravitational softenings 
are taken to be 1.5 kpc for gas and star particles, and 3 kpc 
for dark matter. The total number of particles typically 
becomes ~4x 10 4 at the end of the simulations (z = 0). 

2.3. Data Reduction 

In order to compare the properties of observed and sim- 
ulated galaxies, we compute the observables such as lu- 
minosity and the line-width for each "spiral galaxies" at 
z = 0. Stellar luminosity is computed with the simple 
stellar population (SSP) synthesis models of Kodama & 
Arimoto (1997). The models provide integrated spectra 
along stellar isochrones corresponding to appropriate age 
and metallicity. We take the Salpeter's IMF. Total mag- 
nitude Mi of simulated galaxies is computed to sum up 
star particles using appropriate SSP tables, which are se- 
lected in accordance with age and metallicity of the star 
particles. The line-width W20 is derived similarly to the 
observable as constructing a line-profile of gas weighted 



by mass, and measuring the width at 20% level of a peak 
flux. The catalog of all the simulated galaxies is presented 
in Table §. 

3. RESULTS 

3.1. Evolution and Structures of Individual Galaxies 

Figure [l] shows snapshots of star particles at four red- 
shifts for the case of M = 4 x 10 n M o and A = 0.06. These 
view angles provide a face-on projection at z = 0. Two 
clumps at z = 4 are merging and form a bulge-like system 
between z = 4 and 3. After the bulge formation, the sur- 
rounding gas gradually cools and falls to form a gaseous 
disk, and then a stellar disk is gradually formed. The 
panel at z = 2 shows that a disk-like structure surrounds 
the central bulge. Since the gas density of the inner disk 
is higher than that of the outer part, the star formation 
timcscale is shorter in the inner part than in the outer 
part. The stellar disk is gradually formed from the inside 
to outside in the gas disk. The panel at z — 1 shows a 
larger disk than that of z — 2. 

In Figure || we plot the total star formation rate, total 
magnitude and color in a rest frame against the look-back 
time (and also the redshift). The star formation rate peaks 
at z ~ 3 at a rate of 29M yr _1 , and then declines and 
reaches an almost constant value of lMoyr -1 at z ~ 1, 
when the stellar mass has already been about 80% of the 
total disk mass. At z = 0, about 10% of the total baryonic 
mass remains in the gas and the rest is already in stars, 
which is consistent with the observations of spiral galaxies. 
Total magnitude also peaks near the time of the maximum 
star formation. Then the magnitudes in all passbands are 
monotonically declining until z = because the mean age 
of the stellar component becomes older and massive stars 
gradually die. The amount of the decline between z = 3 
to the present is 1.6 mag in /-band. 

We show the final snapshots of star particles in Fig- 
ure |^. The /-band luminosity profile and rotation- velocity 
profile of the final galaxy are shown in Figure ||. Observed 
internal structures, such as the exponential light-profile 
and flat rotation curve, are well reproduced as the results 
of the simulation. We confirm that all the 14 simulated 
"spiral galaxies" attain such properties. 

About 30% of total gas angular momentum transfers to 
that of dark matter in the period z — 5 — 2, when con- 
densed gas core with dark halos merge into a larger object 
(Katz & Gunn 1991; Navarro, Frenk & White 1995). Af- 
ter the period, the baryon disk scarcely looses the angular 
momentum, because of no infall of lumpy clumps in our 
simulations. The specific angular momenta of the baryon 
disks lie in the range of those in observed spiral galaxies 
at z = (Fall 1983; Contardo et al 1998). In the following 
discussion, we investigate a statistical property of 14 simu- 
lated galaxies, all of which have internal structures similar 
to observed spiral galaxies as discussed above. 

3.2. The Tully-Fisher Relation 

In Figure ^, we compare the observed (open squares) 
and simulated (closed circles) TF plot in /-band. We use 
the observed data presented by Han (1992) and translate 
them to absolute values with distances assuming h = 0.5. 
The slopes of the solid and dashed lines are derived by 
fitting to observed data (Giovanelli et al. 1997a), and the 
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zero points are fitted by eye. In Figure o, we notice that 
the simulated TF relation shows the following three points: 
(i) The slope of TF is well reproduced, (ii) The scatter of 
TF is also similar to the observations, (iii) The zero point 
is systematically fainter. These points are consistent with 
the previous results which included environmental effects, 
e.g., tidal field and infall/outflow of mass (Steinmetz & 
Navarro 1999). Here, we artificially control initial condi- 
tions of isolated spheres, and examine how final structures 
of "spiral" galaxies depend on the initial conditions. Based 
on these simulations, we discuss the TF relation in detail 
on three points: the zero point, slope and scatter. 

The zero point: Our simulated galaxies are 1.5 mag 
fainter than observed ones. We could find a possible so- 
lution of this discrepancy by changing the adopted cos- 
mological model: Changing h (fio) would cause a verti- 
cal (horizontal) shift of the plotted points in Figure |[ 
If larger h (> 0.5) is taken, distance estimations for the 
observed galaxies become smaller, and open squares (ob- 
served galaxies) vertically shift downward in Figure |[ 
Larger h also provides younger age, and therefore, brighter 
luminosity for galaxies, because galaxies loose luminosity 
with their age. It results in an upward shift of the sim- 
ulated galaxies. If smaller Qq (< 1) is taken, the baryon 
fraction fib / fio becomes larger because fib should be fixed 
by the Big Bang Nucleosynthesis. Since /-band luminos- 
ity is almost determined by the total baryonic mass, it 
results in a smaller mass-to-light ratio M/L, which shifts 
the filled circles (simulated galaxies) leftward in Figure |[ 
If we take a cosmological model with larger ft, > 0.5 and/or 
smaller fio < 1, which is suggested by recent observations 
(Giovanelli et al. 1997b, Perlmutter et al. 1998), the simu- 
lated TF relation would become consistent of the observed 
one. Note that the zero point also depends on the redshift 
when the whole mass is merging into a galaxy in hierar- 
chical cosmology (Mo, Mao & White 1998). The fi < 1 
universe will shift the zero point brighter (upward), and 
also alleviate the discrepancy. 

The slope: As seen in Figure ^, our simulations well 
reproduce the slope and scatter of the TF relation. The 
shifted zero point may not affect the mutual comparison 
among the simulated galaxies. In the left panel of Figure 
|^, the simulated galaxies with different masses are repre- 
sented by different symbols on the TF plot. The slope of 
the solid line is derived from observed data (Giovanelli et 
al. 1997a). Figure || clearly shows that larger total mass 
of galaxies result in brighter luminosity and larger line- 
width. The resultant slope of model galaxies on the TF 
plot is quite similar to the observations. The TF correla- 
tion originates in the differences of total masses of galaxies. 

The scatter: Figure || (left panel) shows that J-band lu- 
minosity are different among the galaxies with the same 
mass. The reason of the non-constant M/L is found in 
the right panel of Figure ||. The right panel is the same 
plot as the left panel, but different symbols are used for 
different initial spin parameters A. It shows that galax- 
ies with different A distribute nearly perpendicular to the 
TF correlation, resulting in the scatter of the TF relation. 
The amplitude of the scatter depends on the range of A. 
We took A = 0.04 — 0.10, which is a range if angular mo- 
menta are caused by the cosmological tidal field (Barnes 
& Efstathiou 1987). The role of spin parameters on the 
scatter can be explained in two ways: (i) Low (high) A 



leads to a centrifugally concentrated disk, which provides 
large (small) line- width W^o- (h) Low (high) A results in 
smaller (larger) disk, and high (low) surface density of the 
gas disk, which controls the time-scale of star formation 
TV For low A, is relatively short due to the high surface 
density, and most of stars are formed at the earlier phase 
of galaxy formation. Therefore, the mean age of the stellar 
component in low A galaxies is older than that in high A 
galaxies at z = 0. Low A makes fainter galaxies than high 
A through the age of their stellar component. This effect 
also makes the scatter of the TF relation. 

4. SUMMARY AND DISCUSSION 

We have investigated the origin of the TF relation by 
simulating formation of 14 galaxies with the iV-body/SPH 
direct calculations, which includes the cooling, star for- 
mation and stellar feedback. Internal structures of spiral 
galaxies, such as the exponential disk profile and flat rota- 
tion, are well reproduced as the results of the simulations, 
while these properties were assumed a priori in past semi- 
analytic approaches (Mo, Mao & White 1998; Heavens & 
Jimenez 1999). We observed the total magnitude and gas 
line- widths of the simulated "spiral galaxies", and found 
that they also reproduce the slope and scatter of the TF 
relation, except for the zero-point. We have simulated 
galaxy formation from many different initial conditions, 
and found that the slope of the TF relation is produced 
by the difference of total galactic masses, and the scatter 
is produced by the difference of initial spin parameters. 

We should comment upon our model of star formation 
and feedback. The following issues about star formation 
and feedback remain for future investigations with higher 
resolution simulations. Luminosity in red passbands is de- 
termined firstly by the stellar amount, and secondly by 
their age. Simulated galaxies have the stellar amount sim- 
ilar to the observations, typically ~ 80 — 90% of the to- 
tal baryonic mass at z = 0. Their luminosity hence de- 
pends on stellar age distribution, i.e., star formation his- 
tory. Our star formation model is based on the Schmidt's 
law (p* cx /Og a g), which is an optimal but robust model if 
the mean rate of star formation can be determined by a 
mean local density of the gas. This averaging treatment 
of star formation would be most suitable for the current 
resolution ~ 1 kpc. In the Milky Way Galaxy, star forma- 
tion is observed at molecular clouds, and further discus- 
sion about different kinds of star formation model would 
requires to resolve at least the scale of molecular clouds 
<; lOOpc, 10 times smaller than our spatial resolution. In 
our simulations, the star formation timescale is mostly set 
to 20rd yn . Since the local dynamical timescale Td yn is much 
shorter than galaxy age, the gas is transformed to stars ef- 
ficiently in high redshift. Taking the model which gives 
longer timescales of star formation would provide brighter 
galaxies because galaxies which consist of younger stars 
are brighter. 

Stellar feedback would control the self-regulation mech- 
anisms in spiral galaxies as discussed by Silk (1997). Our 
adopted model behaved like a minimal feedback, because 
feedback energy released as thermal energy which could 
be radiated away in a short timescale. If the energy is re- 
leased as kinetic form as discussed in Navarro & White 
(1993), the feedback would become more efficient and 
could change the history of star formation. It however 
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is not clear that feedback through kinetic energy can be 
occurred in the scale of ~ lkpc, i.e., our spatial resolution. 
Since supernovae remnants expand to about 10 — lOOpc in 
the inter stellar medium, 10 — 100 times higher spatial res- 
olution than the present one should be required for the 
detailed modeling of the local stellar feedback. These is- 
sues remain for future inves tigat ions. 

As discussed in Section 3.2, the discrepancy of the 



1.5 mag fainter zero point in J-band would be alleviated in 
the cosmological model of h > 0.5 and fi < 1. We should 
note, however, that this discrepancy cannot be solved in 
the h = 0.5 and = 1 cosmology with different models 
of star formation and feedback. The /-band luminosity 
declines by only 1.6 mag from z — 3, the time of the peak 
luminosity, to the present. In order to shift the zero point 
1.5 mag brighter by changing the star formation history, 
almost all the stars must be formed at very low redshift 
z < 0.1 (recent 2 Gyr) in all the spiral galaxies in the local 
universe. 

Our initial conditions consider collapse of isolated 
spheres. Under these idealized initial conditions, the infall 
of sub-galactic clumps at z <^ 2, which is expected in the 



fi = 1 CDM universe, cannot occur, though a continu- 
ous gas supply would enhance the star formation rate at a 
low redshift. If we consider later gas infalls, the resultant 
galaxies would have brighte r lu minosity. Note however 
that as discussed in Section 3.2, our derived TF relation 



was consistent with the previous simulations (Steinmetz & 
Navarro 1999), which included the environmental effects, 
e.g., tidal field and infall/outflow of mass. In our isolated 
conditions, the simulations do not suffer from an extreme 
transfer of angular momentum from baryon to dark mat- 
ter (Navarro, Frenk & White 1995; Navarro & Steinmetz 
1997) because no merging occurs at a low redshift, 
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Fig. 1. — Snapshots of star particles at various rcdshifts in the case for M = 4 X 10 11 Af Q and A = 0.06. The origin of each figure is shifted 
to make the object nearly at the center. The view angles gives the face-on projection at z = 0. Note that these plots do not actually show 
the mass density because each star particle does not necessary have the same mass. 
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Time (Gyr) 

Fig. 2. — Evolution of total star formation rate, total magnitude and color in the case for M = 4 X 10 11 Mq and A = 0.06. Upper panel: 
star formation rate vs time. Middle panel: U, B, V, R and I-band total magnitude in a rest frame vs time. Lower panel: V-I, B-U and U-V 
color in a rest frame vs time. The scale in rcdshift is also illustrated in the upper panel. 
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Fig. 3. — Snapshots at z 
panel: a edge-on projection, 
the same mass. 



X (kpc) 

= of star particles in the case for M = 4 X 10 Mq and A = 0.06. Upper panel: a face-on projection. Lower 
Note that this plot does not actually show the mass density because each star particle does not necessary have 



8 



Origin of the Tully-Fisher Relation 




5 10 15 20 

R (kpc) 

Fig. 4. — Internal structures of final galaxies at z = in the case for M = 4 X lO^M© and A = 0.06. Upper panel: the /-band surface 
brightness p-i profile. Lower panel: the rotation curve. All the simulated galaxies finally show the exponential light-profiles and flat rotation 
curves, similar to observed spiral galaxies. 
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Fig. 5. — The TF relation of observed (open squares) and simulated galaxies (filled circles). The slopes of the solid and dashed lines are 
derived by fitting to observed data, and the zero-point are fitted by eye. Simulated galaxies well reproduce the slope and scatter of the TF 
relation except for their 1.5 mag fainter zero point. 
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Mass Angular momentum 




log(W 20 (km/s)) log(W 20 (km/s)) 

Fig. 6. — The TF relations of simulated galaxies. The slopes of the solid lines are derived by fitting to observed data, and the zero-point is 
fitted by eye. Left panel: mass dependence of galaxy distribution in the TF plot. Different masses are explicitly shown with different symbols. 
The TF relation is extended to the direction which traces different masses. Right panel: angular momentum dependence of galaxy distribution 
in the TF plot. Different spin parameters A are shown with different symbols. Galaxies with different A distribute nearly perpendicular to the 
TF relation, and result in the scatter of the TF relation. The amplitude of the scatter depends on the range of A. The range A = 0.04 — O.fO, 
provided by the cosmological tidal field, results in the appropriate scatter. 
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M rls 


Mz r \ s Period 




(erg/star) 


(M /star) 


(M /star) (10 8 yr) 


sw 


10 49 




0-0.4 


SNIa 


10 51 


1.4 


1.4 0-0.4 


SNII 


10 51 


15.3 


2.5 5-30 



Table 1 

A STAR WITH > 8Mq IS ASSUMED TO RELEASE THE TOTAL ENERGY E T i a , MASS M ris AND METAL M Z lls PAR A STAR THROUGH 
STELLAR WIND (SW), TYPE IA (SNIa) AND II (SNII) SUPARNOVAE. THE VALUES LISTED ABOVE ARE AVERAGED OVER STARS 
DISTRIBUTED ACCORDING TO THE SALPETER'S IMF, USING THE YIELD TABLES OF NUCLEOSYNTHESIS (NOMOTO ET AL. 1997A & 

1997b). We set the upper (lower) mass cut of the IMF to 60M q (O.1M ). The number of stars with > 8M in a 

star particle can be estimated with the Salpeter's IMF. 



Run 


M 

(M ) 


A 


Sp/p 


W 20 
(kms -1 ) 


Mi 
(mag) 


1 


8 x 10 11 


0.10 


1.9(7 


466 


-22.31 


2 




0.08 


1.9a 


524 


-22.07 


3 




0.06 


1.9a 


549 


-21.90 


4 


4 x 10 11 


0.08 


1.7a 


420 


-21.32 


5 




0.06 


1.7a 


437 


-21.21 


6 




0.04 


1.7a 


454 


-21.06 


7 


2 x 10 11 


0.10 


1.5a 


299 


-20.78 


8 




0.08 


1.5a 


331 


-20.66 


9 




0.06 


1.5a 


342 


-20.47 


10 




0.04 


1.5a 


349 


-20.34 


11 


1 x 10 11 


0.10 


1.4(7 


236 


-20.09 


12 




0.08 


1.4(7 


248 


-19.92 


13 




0.06 


1.4(7 


257 


-19.80 


14 




0.04 


1.4(7 


254 


-19.66 


Table 2 



A catalog of the simulated galaxies. Our two initial free parameters, total mass M and spin parameter A, 

INITIAL OVERDENSITY 8p/p ABOVE THE BACKGROUND FIELD AND TWO OBSERVABLE VALUES OF FINAL GALAXIES, LINE-WIDTH 

W20 AND /-BAND TOTAL MAGNITUDE Mi, ARE LISTED. 



